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ABSTRACT 

Magnetization in highly conductive plasmas is ubiquitous to astronomical systems. Flows in such 
media can be described by three path functions A a , or, for a steady flow, by two stream functions 
A„ and an additional field such as the mass density p , velocity v, or travel time At. While typical 
analyses of a frozen magnetic field B are problem-specific and involve nonlocal gradients of the fluid 
element position x(t), we derive the general, local (in A or A space) solution B = (dx/dA a ) t B a p/p, 
where Lagrangian constants denoted by a tilde are directly fixed at a boundary hypersurface H 
on which B is known. For a steady flow, pB/p = (dx/dX^A t B K + vB^/v; here the electric field 
E ~ (t? 2 VAi — i?iVA 2 )/p depends only on A K and the boundary conditions. Illustrative special 
cases include compressible axisymmetric flows and incompressible flows around a sphere, showing 
that viscosity and compressibility enhance the magnetization and lead to thicker boundary layers. 
Our method is especially useful for directly computing electric fields, and for addressing upstream 
magnetic fields that vary in spacetime. We thus estimate the electric fields above heliospheres and 
magnetospheres, compute the draping of magnetic substructure around a planetary body, and demon¬ 
strate the resulting inverse polarity reversal layer. Our analysis can be immediately incorporated into 
existing hydrodynamic codes that are based on stream or path functions, in order to passively evolve 
the electromagnetic fields in a simulated flow. Furthermore, in such a prescription, the electromagnetic 
fields are frozen onto the grid, so it may be developed into a fully magnetohydrodynamic (MHD), 
efficient simulation. 

Subject headings: magnetohydrodynamics (MHD) - magnetic fields - planets and satellites: magnetic 
fields - galaxies: magnetic fields - ISM: magnetic fields. 


1. INTRODUCTION 

Magnetic fields frozen in highly conductive plasmas are 
studied in diverse astronomical systems. Examples in¬ 
clude the solar corona, flares and wind (e.g., Longcope 
2005; Zhang & Low 2005), planetary magnetospheres 
and bow shocks (e.g., Spreiter & Alksne 1970; Spreiter 
& Stahara 1995; Zhang et al. 2004; Corona-Romero & 
Gonzalez-Esparza 2013), the interstellar medium (ISM) 
(e.g., Price & Bate 2008; Li et al. 2011; Padovani et al. 
2014), in particular where it meets the solar wind (Parker 
1961; Aleksashov et al. 2000; Whang 2010), and the in- 
tergalactic medium (IGM) of galaxy groups and galaxy 
clusters (Bernikov & Semenov 1979; Kim et al. 1991; 
Vikhlinin et al. 2001; Keshet et al. 2010; Briiggen 2013). 
For general reviews of magnetization in astronomical sys¬ 
tems, see Widrow (2002), Vallee (2011). 

In such systems, magnetization away from shocks and 
reconnection regions can be typically approximated using 
ideal magnetohydrodynamics (MHD), on a background 
ranging from a simple axisymmetric flow around a blunt 
object to complicated, turbulent motions. The passive 
evolution of magnetic fields frozen in a given flow is thus 
important in space physics, astrophysics, applied math¬ 
ematics, and computational physics (e.g., Ranger 1997; 
Sekhar 2003; Sekhar et al. 2005; Bennett 2008). How¬ 
ever, present derivations of the magnetic field evolution 
in a general flow are inherently nonlocal. 

Formally, the ideal MHD equations can be solved to 


give (Elsasser 1956) 

B = ^(B-V)x, ( 1 ) 

where B is the magnetic field, p is the mass density, x 
is the position of the fluid element, and a tilde denotes 
(henceforth) a Lagrangian constant evaluated on a ref¬ 
erence hypersurface H , on which B is known. Thus, 
V is a non-local operator, acting in the vicinity of H, 
rather than around x. Therefore, to derive B from Equa¬ 
tion (1), one must essentially integrate first over the flow 
to compute the mapping of x on H. 

Here we show that B and E can be locally and simply 
computed for a general flow, by working in the space 
spanned by the path functions or stream functions that 
describe the dynamics. Such functions exist as long as 
particle diffusion can be neglected (Yih 1957). Stream 
functions describe steady flows, and are also known as the 
Euler or Clebsch potentials of the mass flux (Euler 1769; 
Clebsch 1859). Path functions describe time-dependent 
flows, and are sometimes referred to as material functions 
(Van Roessel & Hui 1991; Loh & Hui 2000). 

Thus, the dynamics of time dependent flows can be 
fully described by three path functions, A a= i i2j3 (Yih 
1957, in three spatial dimensions; used henceforth and 
denoted 3D). The dynamics of steady flows can be fully 
described by two stream functions, A k= i i2 , which fix the 
mass flux j = pv (Giese 1951), and an additional field tp 
such as p , the velocity v , or the travel time At. (Here we 
ignore an additional thermal quantity, such as temper- 
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ature or pressure, which is needed to fully describe the 
flow but is unnecessary for our purposes.) Path (stream) 
functions are useful in picturing the flow, because their 
equivalue surfaces intersect at pathlines (streamlines), 
the trajectories of fluid elements in spacetime (in space). 

For simplicity, we focus on non-relativistic flows, with 
negligible particle diffusion, in the ideal MHD limit. Ar¬ 
bitrary flows are considered, including those involving 
discontinuity surfaces such as shocks. Our assumptions, 
in particular ideal MHD, may break down across shocks 
and in reconnection regions, where kinetic plasma ef¬ 
fects can be important. In such cases, our analysis may 
be piecewise applied, for example, on both sides of the 
shock. 

The solution we derive for passive magnetization holds 
for arbitrary initial magnetic field configurations. It is 
analyzed for a general steady, axisynnnetric flow, tested 
for an arbitrary time dependent, one-dinrensional (ID) 
flow, and illustrated for basic flow patterns, in particular 
steady flows around blunt objects used to model astro¬ 
nomical systems. This includes both novel solutions and 
generalizations of previously solved problems. The re¬ 
sults indicate that passive magnetization is in general 
enhanced by viscosity and compressibility effects, lead¬ 
ing to stronger magnetization in front of moving objects, 
and thus to thicker boundary layers. 

Analytic solutions to the nonmagnetized flow around 
an obstacle are available only for a handful of time- 
independent cases, such as incompressible, potential 
flows around simple objects ( e.g ., Landau & Lifshitz 
1959; Leal 2007), or the approximate, compressible flow 
in front of an axisymmetric body (Keshet & Naor 2014). 
Solutions to the induced magnetization are even more 
rare, notable examples including the steady, incompress¬ 
ible, potential flow around a sphere (Chacko & Hassanr 
1997; Lyutikov 2006; Dursi & Pfronnner 2008; Romanelli 
et al. 2014) or around simplified surfaces representing 
for example bow shocks (Corona-Romero & Gonzalez- 
Esparza 2013) or the heliopause (Roken et al. 2014). We 
are unaware of a previous magnetization solution even for 
the simple, Stokes (creeping) incompressible steady flow 
around a sphere. Such solutions are derived as special 
cases of our general axisymmetric result. 

Our primary goal is to show that the evolution of the 
electromagnetic field in a given system becomes trivial 
(transparent) if the time-dependent (steady) flow is pa¬ 
rameterized in terms of path (stream) functions. This en¬ 
ables us to quickly reproduce known analytic results, and 
to analytically solve the field evolution in simple flows 
and test problems of astrophysical importance. More¬ 
over, our method can be incorporated into existing nu¬ 
merical simulations that utilize a stream function or path 
function description; this is straightforward for weak 
electromagnetic fields, but appears feasible for full MHD, 
at least in two-dimensions (2D). 

The paper is organized as follows. In Section 2, we de¬ 
rive the temporal evolution of the magnetic field, for both 
steady and time dependent flows. Arbitrary axisymmet¬ 
ric flows are analyzed in Section 3, focusing in particular 
on the axis of symmetry. In Section 4, we present some 
basic flows and analyze their magnetization. In Section 
5, we apply the analysis to a few illustrative astronom¬ 
ical systems. The results are summarized and discussed 
in Section 6, where we also examine the integration of 


our results into simulations that are based on path or 
stream functions, and the prospects for developing a full 
MHD simulation based on this approach. 

2. EVOLVING ELECTROMAGNETIC FIELDS 

We begin in Section 2.1 by introducing the MHD equa¬ 
tions, and switching to a Lagrangian perspective. Next, 
we analyze the evolution of B for arbitrary steady (in 
Section 2.2) and time-dependent (in Section 2.3) flows. 
To demonstrate the consistency of these two frameworks, 
in Section 2.4, an arbitrary steady flow is analyzed using 
the time-dependent formalism. 

2.1. Two Viewpoints of the Relation B oc pi 

Let l be an infinitesimal vector adjoining nearby fluid 
elements that lie on the same magnetic field line (hence¬ 
forth: length element). The magnetic field can be shown 
to evolve in proportion to pi {e.g., Landau & Lifshitz 
1960). We derive this relation in an Eulerian picture 
in Section 2.1.1, and in a Lagrangian framework in Sec¬ 
tion 2.1.2. 


2.1.1. Eulerian picture 
The ideal MHD equations, 


Gauss’ law, 

V B = 0 ; 

(2) 

continuity, 

§Uv-(p.) = 0 ; 

( 3 ) 


dB 


convection, 

— = V x {v x B) ; 

( 4 ) 

and Ohm’s law, 

E = —— x B , 

( 5 ) 


can be combined {e.g., Landau & Lifshitz 1960) to give 
the Helmholtz equation 


d_ 

dt 




( 6 ) 


Here, E is the electric field, and c is the speed of light. 
As l is fixed to the flow, it satisfies the equation 


f t ={ l -V)v. (7) 

Comparing Equations (6) and (7) indicates that as {B/p) 
and l evolve, their ratio remains constant, such that 


B = B?l . (8) 

pi 

This equation holds for both steady and time-dependent 
flows, and even across surfaces of discontinuity. Using 
Equation (8) to evolve B from H guarantees that the 
MHD equations are satisfied everywhere, assuming con¬ 
tinuity and that Gauss’ law holds on H\ see appendix 
A. 

In our non-relativistic, ideal MHD limit, E is small, 
and may be estimated from Equation (5) once B is 
known. 


2.1.2. Lagrangian picture 

For the subsequent analysis, it is useful to visualize the 
magnetic field amplitude as proportional to the density 
of Held lines, i.e. inversely proportional to the distance if 
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between them. Therefore, B in the direction l is inversely 
proportional to the perpendicular (to l) area element A 
spanned by two such distances, ij)\ and xl> 2 . Hence, 

B = B^l , (9) 

where A = l- (i/>i x^) is assumed infinitesimally small. 
As the mass pAl of the fluid element in the parallelogram 
spanned by {Z,■*/>i,i/> 2 } is constant, A/A = pi/(pi), and 
we recover Equation ( 8 ). 


2 . 2 . Steady Flow Magnetization 
2.2.1. General analysis 

The flow can be fully determined by specifying two 
stream functions, A K where k € {1,2} (henceforth), and 
an additional scalar field <p (Giese 1951). Two surfaces, 
defined by constant values of the stream functions, inter¬ 
sect along a streamline, so 

f • VA K = 0 . (10) 


This equation, along with the time-independent continu¬ 
ity equation 

V-(pt7)=0, (11) 

are identically satisfied by the A gauge (Yih 1957) 


x\ x 2 x 3 


hz = (VAi x VA 2 ) 
pv 


d\i 

d\i 

<9At 

dxi 

8 X2 

9X3 

8 X 2 

dX 2 

9X2 

dxi 

9X2 

9X3 


( 12 ) 


where the determinant, as written without Lame coef¬ 
ficients, is valid only in Cartesian coordinates. We use 
stream functions (and later, path functions) as coordi¬ 
nates. For convenience, we take them in units of length 
by introducing the overall normalizations denoted by 
breve (henceforth; e.g., p/p is a dimensionless mass den¬ 
sity). Equation (12), which is known as the Clebsch rep¬ 
resentation for the mass flux (e.g., Stern 1967), indicates 
that pv cannot be chosen as the independent p, and that 
swapping Ai and A 2 would reverse the flow. 

In order to parameterize the position of fluid elements 
in the volume spanned by the flow, we must introduce, 
in addition to Ai and A 2 , a third variable q ( x ), not nec¬ 
essarily equal to p. The surfaces defined by Ai, A 2 and 
q may never overlap; namely, the vectors 

A 1 = VAr; A 2 = VA 2 ; A 3 = ^ Vq (13) 


cannot be coplanar anywhere in the flow (except at stag¬ 
nation points, where the analysis breaks down; see be¬ 
low). Here, v q = dq/dt = v ■ X?q is the rate of change 
in q along the flow (in units of q/t). The normalization 
factor ( v/v q ) is included to keep the {Aj} basis dimen¬ 
sionless. The freedom in the choice of q will be utilized 
in Section 2.2.2. 

It is useful to introduce the reciprocal basis a q = AJ , 
defined (e.g., Kishan 2007) by a, ■ Aj = Sij, where Sij 


is Kronecker’s delta. Namely, 






(14) 


Here, A = {Ai, A 2 }, and it is understood (henceforth) 
that derivatives with respect to Ai (A 2 ) are taken at a 
fixed A 2 (Ai). Like {A,}, the reciprocal basis too is di¬ 
mensionless, and is nowhere coplanar. Note that the aj 
are not necessarily perpendicular to one another, and 
that the ^-dependent a 12 are not necessarily unit vec¬ 
tors. 

An advantage of the {aj} basis is that the n = {1, 2} 
components of any vector V in this basis are independent 
of q, as 





V ■ Ai 


V ■ VAil 

V = V K a K + V 3 a 3 = 

h 2 

V 3 

— 

v-a 2 

v-a 3 

— 

1 

_ 1 


(15) 

Moreover, for l, the n components l K are conserved along 
streamlines, as we show below. We use (henceforth) 
square brackets to denote vectors in the reciprocal basis, 
which in general is neither orthogonal nor normalized, 
and reserve round brackets for vectors in orthonormal 
(in particular, cylindrical, except in Sections 4.3 and 5.1) 
bases. 

Now, consider two infinitesimally close fluid elements, 
which at some instant are at the locations x(X K ,q) and 

x + l = x(X K + l K ,q + S q ) (16) 

= x (X K ,q) + l K (—+ ) x + °( /2 ) ’ 

where we generalized the Einstein summation rule also 
for stream function indices n = ( 1 , 2 } (and, later on, for 
path function indices a = {1,2,3}). As the A of each of 
the two fluid elements is conserved along the flow, the 
intervals l K = l K are constant along the streamline. 

Henceforth, we work to first order in the infinitesimal 
l, such that 


l — l K a K + l 3 a 3 — l K a K + l 3 a 3 — 


where l 3 = vS q /v q . As l K = l K , the k components of 
l in the reciprocal basis are indeed seen to be constant 
along the streamline. In contrast, l 3 is in general not 
conserved. The evolution of l and its decomposition in 
the reciprocal basis are illustrated in Figure 1. 

In order to evaluate l 3 , consider the travel time At of 
a fluid element from H to the position x , 

q(x) 

At (x) = [ — . (18) 

J V q 

q(x) 


(17) 


During the same time, the x + l fluid element travels 
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from q + S q (not necessarily on H) to q + 6 q , so 


q q+< 5g 

f = At = f _^_ 

l MW) 1 v q (x K +I K ,q') 

<7 q+Oq \ / 


Expanding to first order in l indicates that 



(19) 


( 20 ) 


where derivatives of At = At(X,q,q) with respect to Ai 
or A 2 are taken (henceforth) with the other three argu¬ 
ments of At fixed. The component l 3 can now be related 
to the Lagrangian constants through 


S q j v ~ dAt 
h = v— = < 3 - - vl K —— 
v„ v a X K 


( 21 ) 


If a surface of constant q overlaps with one of the constant 
X K surfaces, v q vanishes and Equation (20) is rendered 
invalid; a piecewise analysis may still be possible. The 
analysis breaks down at stagnation points, where At (and 
so also l and B) diverges. 

Finally, the length element l is now given, throughout 
the flow connected to H , by 


V ~ _dAt\ 

r h ] 


o' 




12 

— 

0 

v V dX.J 

V J 
r? 3 


V 


( 22 ) 

Note that in addition to a K , the At derivative also de¬ 
pends on q. It can be written as 


dAt 

d\: 



(23) 


The Lagrangian constants are fixed at x, where the 
streamline meets H. The point x(x) can be computed 
by inverting H(x), as X(x) = X(x) is given, and re¬ 
mains constant along the streamline. The components 
lj are found by projecting l onto the reciprocal basis as 
in Equation (15), 


1=1- VX K = (l ■ VX k )h and l 3 = (- l ■ Vq) . 

\ V q J H 
(24) 

Equations ( 8 ) and (22) yield the magnetic field, as 


-B = bL 


Bi 


’o' 

b 2 

- 

0 

ib 3 

L v °J 


V 


B, 


dAt 

~dK, 


(25) 


where the Lagrangian constants Bj = ( B/l)lj are fixed 
by the projection Equation (15) on H , 


B k = {B-VX k )„ and B 3 = (— B-Vq] .(26) 

\ v q J h 

Next, consider the electric field E. The vector prod¬ 
uct (v x B) in Equation (5) eliminates the a 3 terms in 



Fig. 1.— Illustration of magnetic field evolution in the flow 
around a blunt object (green half disk), starting from a bound¬ 
ary hypersurface H (dashed cyan curve). A magnetic field line B 
(dotted-dashed red curve) is locally represented by an infinitesi¬ 
mal vector l (thick red arrow), adjoining two fluid elements (pink 
disks) that are initially along B, and advected along streamlines 
(thin black arrows with stream function labels). In the reciprocal 
basis, l is decomposed (see labels) into components parallel to the 
streamline (^ 303 , purple arrows) and along the (constant q) stream 
function gradients blue arrows; here initially confined to H 

for the specific choice q = At; see Section 2.2.2). 


Equation (25), leaving 

£; = ^( j B 2 VA 1 -B 1 VA 2 ) . (27) 

Therefore, E is entirely determined by the stream func¬ 
tions, independent of <p, once B and p are given on H. 
Notice that deriving the electric field directly from the 
stream (and below, from the path) functions is computa¬ 
tionally much simpler than first computing the magnetic 
field, which typically involves integration, and then using 
Ohm’s law to find E. 


2.2.2. Choice of a third variable q(x) 

There are various choices for the third variable, q , de¬ 
pending on the specific problem at hand. A useful choice 
is the travel time At, which in some cases is easily esti¬ 
mated or computed. For q = At, which is the specific 
case illustrated in Figure 1, Equations (22) and (25) re¬ 
duce to 


1 = 1 . 



At 



h 


and 


B = 




B 3 


uj p 


A 

b 2 


p 

p 


(28) 


(29) 


Notice that Equation (28) could have been directly in¬ 
ferred from the reciprocal basis projection Equation (15), 
by noting that l ■ V(A t) = At(x + l) — At(x) + 0(l 2 ) 
is constant (to first order) along a stream line, as this 
is the difference between the travel times of two nearby 
fluid elements. 

In the case of an irrotational flow, a potential $ can be 
defined such that v = where the normalization fac¬ 

tor assigns $ with units of length. Here we may choose 
q = 4>, and read l and B off Equations (22) and (25), 
respectively, with the time derivative term (c/. Equa- 
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tion (23)) 

dAt 

d\ K 



(30) 


Equipotential surfaces are perpendicular to the stream¬ 
lines, so here a\ and a 2 are conveniently perpendicular 
to < 23 . In simple, e.g., axisymmetric, flows where ai and 
a 2 are perpendicular to each other, the reciprocal basis 
thus becomes orthogonal. 

Among the infinite other possible choices of q , notewor¬ 
thy are the natural coordinates of the flow, for example 
the axial coordinate z in simple axisymmetric flows, and 
the length s(x) of the streamline from H to x. In gen¬ 
eral, the simplest variable which is monotonic along the 
flow is advantageous. 


2.3. Time-dependent Flow Magnetization 

A time-dependent flow can be fully specified by three 
path functions, A a= i j 2 i 3 , when particle diffusion can be 
neglected (Yih 1957). Three spacetime manifolds, de¬ 
fined by some constant values of these path functions, 
intersect along a pathline. Therefore, along the flow 

dA a (x,t) ( d 

= (s + ” v ) A " (a! -' > = 0 ' (31) 


Define l as the distance between two simultaneous, in¬ 
finitesimally close spacetime events, 


(A a ,Q) 



(36) 


and 



+ dx^ = (A a + L a , Q + S Q ) (37) 


— x^ T- L a 


dx 11 

dA a 


+ <5 


dx M 

Q W + 


o(i 2 ) • 


As x M = x tJ ’(A a ,Q), derivatives with respect to Q or 
one of the A Q are understood as taken with the other 
three arguments fixed. Here, the infinitesimal L a = L a 
are constant along the pathline, because they are the 
differences between the conserved path function values. 
In contrast, S Q varies in general along the pathline. 

Working to first order in l, the spacetime length ele¬ 
ment is now given by 



f ( dt/dA a 
La ^ dx/dA a 



(38) 


where we defined v Q = dQ/dt as the rate of change in Q 
along the flow (in units of Q/t). The vanishing temporal 
component requires 


which determines v through 

_ € ijk (dtAj) (VAj x VA fc ) . 

2 (VAi x VA 2 ) • VA 3 ’ 1 J 

where is the Levi-Civita symbol. The continuity 
equation (3) is then identically satisfied by the A gauge 
(Yih 1957) 

£ = - (VAi x VA 2 ) • VA 3 . (33) 

P 


S Ji 

Vq 


= -L n 


dt 

dTa ’ 


so the spatial component gives 


l — L a 


dx dt \ 

dA^~ V dKj 


(39) 


(40) 


Taking the scalar product of this result with VA Q , and 
using Equation (31) and the chain rule, implies that the 
constant L a satisfy 


Equations (32) and (33) can be written as (Yih 1957) 


P 

P 



(1 + «) 


Xi 

X 2 

X 3 

1 

dAi 

8 A i 

9Ai 

dAi 

dxi 

dx 2 

dx 3 

dt 

dA 2 

dA 2 

9A 2 

dA 2 

dxi 

dx 2 

8 x 3 

dt 

d \ 3 

dA 3 

0A 3 

dA 3 

dxi 

dx 2 

dx 3 

dt 


where 1 is a unit vector in some non-spatial dimension, 
and the determinant, as written without Lame coeffi¬ 
cients, holds only in Cartesian coordinates. Equivalently 
(Van Roessel & Hui 1991), if v M = (l,u), 


P a = _ apis dAidA 2 dA 3 
p l 6 dx^dx^dx 6 ’ 


(35) 


where {a, /3, 7 ,5} are four-indices, a ; 0 = t, and is 

the four-dimensional Levi-Civita symbol. 

To parameterize spacetime using the path functions 
Ai, A 2 , and A 3 , we normalize them in units of length, 
and introduce a fourth coordinate, Q(x,t). The fourth 
coordinate can be chosen arbitrarily, as long as the equiv- 
alue hypersurfaces of the four parameters never overlap. 


L a =l ■ VA a = (l ■ VA q )^ , (41) 

and are independent of Q. They may thus be determined 
at any position where l is known, in particular on H. 
Finally, combining Equations ( 8 ) and (40) yields 


B — B a 


P 

P 


dx dt \ 
dA a dA a ) 


(42) 


where the Lagrangian constants B a = ( B/l)L a can be 
determined using 


B a = - p B- VA a = (B ■ VA a ) k . (43) 

As in the time-independent case, here too there is con¬ 
siderable freedom in choosing Q. A useful choice is the 
temporal coordinate, Q = t, for which Equation (40) 
reduces to 


l — L a 



(44) 


The time-dependent vectors (dx/dA a ) t can now be iden¬ 
tified as the generalizations of the reciprocal basis vec¬ 
tors found for steady flows; cf. Equations (14) and (17). 
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Equation (42) now gives the magnetic field, 


B = B n 



(45) 


Other choices of Q may be advantageous in certain cir¬ 
cumstances. 

The electric field in a time-dependent flow becomes, 
according to Equations (5) and (42), 

£ = e «‘fKi5r) VA ‘- (46) 


2.4. Time-dependent Description of a Steady Flow 

To ascertain that Section 2.2 and Section 2.3 are mu¬ 
tually consistent, we consider an arbitrary steady flow, 
and regard it as time-dependent. If the steady flow is 
characterized by the stream functions Ai and A 2 , we can 
choose two path functions as time independent, 

Ai = Ai and A 2 = A 2 . (47) 


Solving Equations (32) and (33) with the aid of Equa¬ 
tion (12) and the definitions of Section 2.2, we find that 
the third path function must be 

A 3 = v (t — At) . (48) 

We can now use the time-dependent results of Sec¬ 
tion 2.3, and attempt to reproduce the steady flow results 
of Section 2.2. 

For simplicity, we choose Q = t, such that the solution 
Equation (44) for l in the time dependent regime gives 


i = l k 



+ L 3 



A «,t 



(49) 


This is formally identical to Equation (28), derived for a 
steady flow under the choice q = At, provided that one 
may consistently identify l K = L K and l 3 = —L 3 (v/v). 
Indeed, Equation (41) yields L K = / • V A K = l- VA K = l K , 
and, on H, L 3 = l VA3 = -vl-'V(At) = ~(v/v)(v/v q )l- 
Vq = —{v/v)h; consistent with the definitions of Equa¬ 
tion (24). 

In conclusion, the time-dependent analysis with Q = t 
yields the time-independent analysis with q = At, when 
applied to a steady flow. This shows that the two for¬ 
malisms, involving path functions and stream functions, 
are consistent with each other. 

With the above A a , the time-dependent electric field in 
Equation (46) directly reduces to the time-independent 
E in Equation (27), as dAj/dt = vS 3 j. 


3. STEADY AXISYMMETRIC FLOW MAGNETIZATION 
3.1. General Analysis in Cylindrical Coordinates 

Consider a steady, axisynnnetric flow. Using cylindri¬ 
cal coordinates {p, <f>, z}, such flows are characterized by 
Vif, = 0, and by vanishing f derivatives of the flow param¬ 
eters. Therefore, constant <f> surfaces are perpendicular 
to the flow, and we may choose Ai = Rf, where R is some 
arbitrary length scale of the problem. Since the velocity 


field v is ^-independent, we may choose A 2 = A (g,z). 
For these stream functions, Equation (12) yields 

t. = K,0,»,) = i(g,0,-g) , (50) 

where fi = pgv is the mass flux in a ring of radius g 
(and has units of viscosity). Here, the normalizations R 
and fi = pv R arise from our use of stream functions as 
coordinates with units of length. Round bracket vectors 
in this section pertain to cylindrical coordinates. 

In addition to Ai and A 2 , we must now choose a third 
variable q to parameterize space. The choice q = z is 
advantageous for simple flow regions in which v z does 
not change sign. For our choice of stream functions and 
q, the reciprocal basis Equation (14) becomes 


Q 2 H ~ ~ ,,,, 

ai = a 2 = - Q ; a 3 = v . (51) 

R Hz 

A given B = (B e , B^, B z ) may now be decomposed as 
Bid,i, where (cf. Equation (26)) 



Br' 


T) R 

B<PJ _ 

B = 

b 2 

— 

(. HqB z — h z Bq)/h 


b 3 


1 _ 


(52) 


Recall that in our notations, vectors in square brackets 
(as in Equation (52)) are in the reciprocal basis. 
Equation (25) then yields 


?B = iB^- ^B 2 q+1(b 3 -B 2 d^) (53) 

PR Hz v \ a\ ) 


( &LB p \ 


/ye _ 


(v\ 

fJ-Z S 


Vz Mz 



f B* 

+ 

0 

B z - 

0 



l 0 ) 


{Vz) 


The only computation necessary is an integral for the 
time derivative term, 


dAt 

~d\ 


<9A t(g, z) 
d\ 


4>,Z,Z 


z 



z 


The electric field is found from Equations (27), (50) 
and (52), and given simply by 





/ 

l^z B(p 

E = - 

c 

&b 2 

Qfl Z 

V 

jlc 

Q 

Q 

{fi e B z fi z B e ) 




\ 

—l^gBff) y 


(55) 


3.2. The Axis of Symmetry 

Next, we focus on the vicinity of the axis of symmetry, 
g = 0, along which v e = 0. Taking the incident flow 
in the (— 2 ) direction, we denote the incoming velocity 
as u(p ~ 0 , 2 ) = — v z , and assume it is finite. Then, 
Equation (50) indicates that in the vicinity of the axis, 
A = ( pu/2fi)g 2 + 0(g 3 ). Along a streamline, A is con¬ 
stant,therefore so is pug 2 -, this is a direct expression of 
mass conservation. In the vicinity of the axis, we may 
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thus write g/g = ( pu/pu ) 1 / 2 along a streamline. This 
holds even in the presence of a shock, which conserves 
the normal, in this case axial, mass flux j = pu. 
Equation (53) now becomes 


irrotational (and thus, effectively inviscid) case in Sec¬ 
tion 4.1, and in the viscous limit in Section 4.2. In Sec¬ 
tion 4.3 we study the simple, ID class of time-dependent 
flows. 






B = 

(jit} 

+ c 

0 





\ ^ J 




where 

Z 

r , \ _ jfj-z dAt _ jit r u\dz' 

(z) ~ m ox ~ vi J 

z 


(56) 


(57) 


and we expanded u = uq(z) + u\(z)p + 0(p 2 ). 

For a potential flow, u\ = 0, and so C = 0. Such a 
flow arises, for example, if a homogeneous incident flow 
remains subsonic or mildly supersonic ( e.g Landau & 
Lifshitz 1959). 

It is typically possible to derive (p/p) as a function 
of u, such that B = B(u ) along the axis. For exam¬ 
ple, consider a steady, inviscid flow in a poly tropic gas 
of an adiabatic index 7 , with dynamically insignificant 
magnetic field. Such a flow is governed, in addition to 
continuity, by Euler’s equation, 


(v • V) v = 


VP 

P 


c|Vp 

P 


(58) 


and Bernoulli’s equation, 


u ct ct 

— H-—- = —— = constant . 

2 7 — 1 7 — 1 


(59) 


Here, P is the fluid pressure, c s = ( 7 P/p ) 1 / 2 is the (local) 
speed of sound, and a bar denotes (henceforth) a putative 
stagnation point where v = 0 , whether or not such a 
point lies along a given streamline. 

The axial mass density solution is p 2 oc ( S 2 — Mq) s , 
where S = [ 2/(7 — l)] 1 / 2 , and M 0 = u/c s is the ax¬ 
ial Mach number with respect to the stagnation sound 
speed. Hence, away from shocks, 


p f S 2 - M 2 \ 1/h ~ 1] 
P ~ \S 2 -M 2 ) 


(60) 


can be used in Equation (56) to find B once u is deter¬ 
mined. When a shock is present, the Lagrangian con¬ 
stants here are taken downstream; if H is located up¬ 
stream, the shock compression ratio ( Mq/W) 2 should be 
incorporated. Here, W = [ 2 /( 7 + l )] 1 / 2 and S are the 
weak and strong shock limits of Mo; see Keshet & Naor 
(2014). 


4. MAGNETIZATION OF BASIC FLOWS 

To illustrate the general magnetization solution, de¬ 
rived in general in Section 2 and for axial symmetry in 
Section 3, here we study the magnetic fields evolving in 
specific, basic flows. Passive magnetization in steady, in¬ 
compressible flows around a sphere are analyzed in the 


4.1. Potential Incompressible Steady Flow Around a 

Sphere 

Consider the inviscid, incompressible, steady flow 
around a solid sphere. We assume that the incident flow 
is homogeneous, unidirectional, and subsonic, so the en¬ 
tire flow is irrotational. Here, the mass density p = p is 
constant, the radius of the sphere is taken (henceforth) 
as R = 1, and the constant v = —uz is taken on H far 
upstream, z —> 00 . It is thus natural to choose p = p, 
R = 1 and v = u, so p = pu. 

As V • v = 0 and V x v = 0, the flow satisfies the 
Poisson equation V 2< I> = 0 for the velocity potential $, 
defined by v = ■aVd*. The solution for boundary con¬ 
ditions of a vanishing normal velocity across the sphere 
is 

®=- z ( i+ 7)' <61) 


Here we primarily use cylindrical coordinates, but occa¬ 
sionally allude to spherical coordinates {r, 6, </>}, where 
r 2 = g 2 + z 2 . Thus, 



and v e 


3ugz 
2 r 5 


(62) 


The stream functions can be chosen as Ai = </> and, ac¬ 
cording to Equation (50), 

, (63) 


such that g = (2A) 1 / 2 . 

Solving Euler’s Equation (58) gives the pressure, 


P = P 


pu' 


4r 3 — 5 + 3 (4r 3 — l) cos 2 9 
16 r 6 


(64) 


where cos 29 = (z 2 — g 2 )/r 2 . Note that an incident pres¬ 
sure smaller than (5/8 )pu 2 would lead to a nonphysical 
P < 0 around z = 0, as the plasma would be too cold to 
remain incompressible. 

Consider the evolution of a magnetic field, initially 
given by B = (B g , B^, B z ), for different choices of q. In 
this subsection, round-bracket vectors are in cylindrical 
coordinates. 


4.1.1. The natural choice q = z 

For the choice q = z, which is natural because z rnono- 
tonically decreases for all fluid elements, the magnetic 
field is given by (cf. Equation 53) 

B = B^ej) — B e ^-Q — ^ (b z + B e gu^^\ , (65) 
g gv z u \ o\ J 

where 

dAt _ fdAt\ _ 12 j r' 8 (r ' 2 - 5z' 2 ) dz' 

~8X ~ \ d\ ) 1 (r ' 2 + 2r ' 5 - 3z ' 2 ) 3 ' 

OO 

(66) 
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Fig. 2. — Potential, incompressible, steady flow (arrows with direction and length proportional to v) around a solid sphere (green half 
disk), for a uniform magnetic field far upstream. The constant cf> slice shows the pressure enhancement (gray scale; for P = pu 2 ) and the 
magnetic amplification B/B (colorscale; Green 2011) for field lines parallel (arrows) or initially perpendicular (solid curves) to the flow. 


The integral is taken at a constant A, so r'( A, z') is found 
from Equation (63). 

Figure 2 shows a constant </> slice through the flow, de¬ 
picting the pressure (grayscale) and the magnetic field 
components parallel (colorscaled arrows) or initially per¬ 
pendicular (colorscaled curves) to the velocity. For illus¬ 
trative purposes, in the figure we take P = pu 2 . 

The electric field here is a special case of Equation (55), 
namely 

/ 


E = 


Q_ 

QC 


\ 


v z B tl 

§* u ^e 

\—VgB<b ) 


4.1.2. The problematic choice q = r 


4.1.3. The orthogonal choice q = <f> 

Another optional choice for q is the velocity potential 
$. Here, the reciprocal basis becomes 


a i — q 4 >; a 2 — 


gv 


(v x 0 ) ; 


a 3 = v 


(70) 


These vectors are perpendicular to one another, and 
therefore may be useful as an orthogonal (albeit not or¬ 
thonormal) basis. The magnetic field here is given by 


(67) 

0 1 B<p 


0 

B = 

QB e 

— 

0 


1 

S|>3 

_ 1 


gv 


B, 


dAt 

’~d\ 


(71) 


Another seemingly possible choice is q = r. According 
to Equation (25), the magnetic field may then be written 
as 

B = B^4> + B g ^6-Z(B z + B e gu^), ( 68 ) 

Q QZ U \ CM J 

where 

dAt _ /dAt\ 

~d\ ~ \ dX ) 


(69) 

However, this result, obtained previously (Dursi & 
Pfrommer 2008; Romanelli et al. 2014) by solving the 
MHD equations, is strictly valid only in the half space 
9 < 7 t/2 . At the surface 9 = 7 r/ 2 , the denominator of the 
integrand diverges, and the integration cannot be contin¬ 
ued. This is to be expected, as surfaces of constant A and 
surfaces of constant r overlap at 0 = 7 r/ 2 , asVA and Vr 
become parallel. The problem can be circumvented by 
computing the converging At, rather than its diverging 
derivative. 

The integral in Equation (69) can be solved analyti¬ 
cally as a power series in A; see appendix B. 


—r ,4 dr'/i 


(r ,s — 1) 2 


where 


dAt _ / 

~d\ ~ V 


$ 



— OO 


dAt \ 

z' A + 4r , 8 4 >' 2 + 10^ , 4 r ' 3 (l - 2r' 3 ) 
r~ u z~ 4 \Ar' & <E 2 + 3 z' 4 (1 — 4r ' 3 )] 3 


(72) 


The parameters r'(A, $') and z'(A, $') are found from 
Equations (61) and (63). 


4.2. Stokes Incompressible Steady Flow Around a 
Sphere 

Consider the Stokes, incompressible, steady flow 
around a solid sphere. As in Section 4.1, p = p is con¬ 
stant, the radius of the sphere is taken as R = 1, we 
assume that the incident flow v = —uz , taken on H far 
upstream (z —>■ 00 ), is homogeneous and unidirectional, 
and use a mixture of cylindrical and spherical coordi¬ 
nates. It is again natural to choose p = p, R = 1 and 
v = u, so p = pu. 

Stokes flows are characterized by viscous forces much 
greater than the inertial forces, so the flow is governed 
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by the approximate Navier-Stokes equation 


which is formally equivalent to Equation (65). Here, 


0 = -V 



+ u'V 2 v, 


(73) 


where v is the kinematic viscosity, assumed constant. 
Taking the curl of this equation gives 

V x (V 2 u) = 0 . (74) 

The solution, for boundary conditions of a hard sphere, 
is 


V z 



1 


(r 2 — 3 z 2 ) (r + 1 ) 
4 r 4 


and 


v e = 3m I 1 


QZ_ 

j'Z J 


(75) 


(76) 


Solving Equation (73) gives the pressure distribution, 


P = P + puv 


3z 

2 ^ ' 


(77) 


Here too, the incident plasma cannot be too cold, as an 
incident pressure P < (3/2)pu 2 /Re, where Re = uR/v 
is the incident Reynolds number, would lead to a non¬ 
physical P < 0 around the z < 0 axis. 

The stream functions can be chosen as Ai = <f> and, 
according to Equation (50), 


A = A 2 


r 

2 



(78) 


As in the case of a potential flow, here too g = (2A) 1 / 2 , 
and we proceed similarly with the natural choice q = z. 
The magnetic field is given by 


B = B^4>-B e 
Q 


QU „ V 

- Q- - 

gv z u 


+ B e gu 



(79) 


dAt, 

~d\ 



(80) 



r' 6 _ 5 ^V 2 + r ' 4 (1 + 3^ 2 ) 
(l-^) 3 (r' 2 -3^ + ^) 3 


dz’ 


where r’{\,z') is found from Equation (78). 

Figure 3 shows the pressure and magnetization in such 
a flow, with the same notations used in Figure 2. For 
illustrative purposes, in the figure P = 2 pu 2 /Re. 

Here too, the electric field is a special case of Equa¬ 
tion (55), and is given by Equation (67). 


4.3. One-dimensional, Time Dependent Flow 

As a simple illustration of magnetization in the time- 
dependent regime, consider the case of an arbitrary ID 
flow. Classical examples include the plasma in a cylin¬ 
drical pipe, disturbed by a moving piston, or a ID rar¬ 
efaction wave. If discontinuities are present in the flow, 
the following applies separately to each region bounded 
by surfaces of discontinuity or H. 

Without loss of generality, take the flow in the x di¬ 
rection, such that v = vx and the perpendicular velocity 
components vanish. A simple choice of path functions is 
then 


A 2 = U and A 3 = z , 


(81) 


as both are constant along the flow. Equations (32) and 
(33) therefore constrain Ai, as 


dAi 

dx 


= — — and 
y,z,t P 


0 Ai 

dt 


pv 

P 


The solution is the path function family 

1 r 

Ai (x,t) =—— / p(x',t) dx' + C , 
P J 


(82) 


(83) 


where we used the continuity equation, and C is an in¬ 
significant constant. 
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For a given l = (l x ,l y ,l z ), Equation (41) yields the 
Lagrangian constants 

L\ = — j;l x ; L 2 = l y j L 3 = l z , (84) 

where here and in the rest of Section 4.3, round brack¬ 
ets are in Cartesian coordinates. Equivalently, for B = 
(B x , By, B z ), Equation (43) gives 

Si = - P -B x ; B 2 =B y ; B 3 = B z . (85) 

Finally, Equations (44) and (45) give, respectively 


Consider however the case of a steady flow, with the 
choice q = At. Here, even if Bi n (x , t) varies in space and 
time, we may directly determine B using 

B(x; t) = B(Xi, X 2 , At; t ) = J5(Ai, A 2 ,0 ;t = t — At) 

= X 2 , 0; t) , (88) 

where f is the time at which the fluid element, which 
presently lies at (x,t), was on H. Equation (88) applies 
for any choice of q , if Equation (18) is used to compute 
the travel time. For example, in cylindrical coordinates 
with H at z —> oo, we may compute it by 


l = , ( 86 ) 

and 

B=(b x ,£b v ,£b z > ) . (87) 

Equations (86) and (87) can be seen directly from La¬ 
grangian considerations. Since the gas flows only in the 
x direction, the perpendicular (to the flow) distance be¬ 
tween fluid elements is constant. Perpendicular length 
elements {l y ,l z } and the parallel magnetic field B x , are 
thus conserved. The parallel distance between nearby 
fluid elements is simply l x oc p _1 , so the perpendicular 
fields follow {B y ,B z } oc p. 

5. ASTROPHYSICAL EXAMPLES 





oo 


An arbitrary unperturbed B(x,t) can be Fourier de¬ 
composed into components of wavenumber k and angular 
frequency w. As a simple example, assume a flow in the 
(—z) direction, and consider a Fourier component of the 
form 

B = Bo(cos(wt), 0, cos (kg)) , (90) 

where Bq is a constant, and round-bracket vectors in this 
subsection are in Cartesian coordinates. Notice that no 
spatial component of this B contributes to Gauss’ law, 
and that a t-dependence on H induces a ^-dependence 
in its vicinity, i.e. Equation (90) becomes 


The simplicity of the above MHD formalism allows 
us to easily obtain results, both general and specific, 
for the evolution of electromagnetic fields in a wide 
range of astronomical systems. In particular, the for¬ 
malism yields analytic expressions for the electric field, 
directly incorporates time-dependent flows, and admits 
inhomogeneous initial conditions in both space and time, 
B — B(x;t), whether or not the flow itself is time- 
dependent. 

In Section 5.1, we demonstrate such inhomogeneous 
scenarios by considering substructure in B. In particu¬ 
lar, some observational consequences are outlined for a 
stellar wind behind a planet or a moon. Near the nose 
of a body, the magnetic field is greatly amplified; the re¬ 
sulting magnetic layer, analyzed in Section 5.2, becomes 
dynamically significant and distorts the flow. In Sec¬ 
tion 5.3, we derive analytic expressions for the electric 
field above heliospheres and magnetospheres. Sophisti¬ 
cated astronomical scenarios typically require a numeri¬ 
cal computation; the integration of our method into nu¬ 
merical simulations is discussed in Section 6. 

5.1. Spacetime Dependent B 

Given a flow and an initial magnetic configuration 
B ln on H , the electromagnetic field solution in Equa¬ 
tions (25), (27), (29), (42), (45) or (46), depends only 
on knowing the Lagrangian constant B corresponding 
to every spacetime location. For an initially uniform 
magnetic field, B is everywhere constant, and the result 
is straightforward even for time-dependent flows. For a 
spacetime-dependent initial magnetic field, this mapping 
can become complicated, especially for time-dependent 
flows. 


B = _B 0 (cos[w(f + z/v)\, 0, cos {kg)) (91) 

in the region where the flow is approximately uniform. 

Figure. 4 presents the evolution of the Fourier compo¬ 
nent Equation (90), in an incompressible potential flow 
around a sphere (see Section 4.1). As Figures 4 and 5 
show, substructure in the far upstream magnetic field is 
non-trivially imprinted on the downstream flow. In situ 
measurements of B downward of an object such as a 
planet would thus be sensitive to both the substructure 
and details of the flow. 

As an example, consider measurements of the inverse 
polarity reversal layer (henceforth IPRL) induced by a 
planetary body (Ronranelli et al. 2014). The IPRL is de¬ 
fined as the layer in which the magnetic field component 
in the direction of the initial flow changes its sign, imply¬ 
ing in our case that B z = 0. According to Equation (53), 
for q = z, the IPRL is described for a general flow by the 
simple formula 


B z = B 2 v 


dAt 

ox 





dAt 

~dX 


(92) 


Consider the IPRL for a magnetic field with initial os¬ 
cillatory substructure. For an initial magnetic inclination 
0 O , let 


B = B 0 (sin 0 O 7 0, cos6*o) + x(B 0 sin (9 0 ) cosuit , (93) 


where £ < 1 is a small number; this B trivially satisfies 
Gauss’ law near H. A snapshot of the resulting magnetic 
configuration is shown in Figure 5, where for illustrative 
purposes we adopt the incompressible, steady potential 
flow of Section 4.1. As the figure shows, the far upstream 




Fig. 4.— Evolution of a magnetic Fourier component, in an incompressible potential steady flow around a solid sphere. Notations are as 
in Figure 2; here arrows are either along or opposite to the flow. The magnetic field far upstream of the sphere is given by Equation (90) 
with k = 2, uj = 7r, and t = 0. 


magnetic substructure with wavelength 27 tu/uj is directly 
imprinted on the IPRL. 

Such a model may describe the magnetization not only 
around planets and moons, but also around stellar sys¬ 
tems moving through the ISM, around clumps and bub¬ 
bles moving through the IGM, etc. 


5.2. Magnetic Layers 

Under our assumptions, as the flow approaches the 
stagnation point in front of an object, the magnetic field 
initially perpendicular to the flow diverges, as seen in 
Figures 2-5. The divergence here, around the object, 
and in the wake behind it, may be regarded as due to the 
accumulation of infinitely many field lines, continuously 
stretched around the object, eventually wrapped closely 
around it and stretching to infinity along the wake. 

These nonphysical divergencies signal the breakdown 
of our assumptions in the above regions. In particular, 
the magnetic field here can no longer be regarded as dy¬ 
namically insignificant, and eventually resistivity can no 
longer be neglected (e.g., Chacko & Hassam 1997; Lyu- 
tikov 2006). 

It is useful to consider the extent of the highly magne¬ 
tized layer forming in front of an object in an arbitrary 
axisymmetric flow. Such layers form as the magnetized 
solar wind drapes around a planetary magnetopause, as 
the magnetized ISM drapes around a bubble or a stellar 
heliopause, and as the magnetized IGM drapes around 
a large scale-clump or bubble. The regions surround¬ 
ing such objects could involve plasma depletion layers 
with dynamically significant magnetic fields, and lead to 
plasma instabilities, magnetic reconnection, and turbu¬ 
lence (e.g., Dibraccio et al. 2013; Sulaiman et al. 2014). 
In the following, we neglect the small modifications to 
the flow far from the object due to the presence of a 
highly magnetized region. 

Define the magnetic pressure fraction 


es = /3 P 1 


B 2 /8n 

Ph 


(94) 


where (3 P is the beta of the plasma, and P th is the ther¬ 
mal pressure. The fractional magnetic amplification may 


then be defined as 



Figure 6 shows the amplifications of magnetic pressure 
(. B 2 /B 2 ) and pressure fraction (cb/cb), as a function of 
the distance S = r — 1 in front of the unit sphere, for 
various steady flows: an incompressible potential flow, 
a Stokes incompressible flow, and inviscid compressible 
flows in the subsonic and mildly supersonic regimes. As 
the above flows are approximately (for the mildly su¬ 
personic case), or precisely, irrotational along the axis, 
Equations (56) and (57) become 

B(z) 2 = P -!^B\ + C^B 2 , (96) 

p/u j 2 11 

where B± and B\\ are the axial magnetic field compo¬ 
nents initially perpendicular and parallel to the flow. 
The figure focuses on the diverging, perpendicular com¬ 
ponent. 

As Figure 6 shows, viscosity induces magnetic amplifi¬ 
cation farther from the body. Compressibility has a sim¬ 
ilar, but somewhat weaker, effect on magnetization, but 
the magnetic pressure fraction is not strongly affected. 
The large difference in es between the incompressible 
and compressible inviscid flows arises from the different 
pressure profiles corresponding to the different equations 
of state. 

Close to stagnation, the flow is adiabatic, p and P th ap¬ 
proach constant values, and B becomes primarily trans¬ 
verse, so 

B 2 „2 

cb oc oc B 2 , oc — oc « _1 , (97) 

Pth J 

which scales as 5~ 2 for the Stokes flow, and as d -1 for the 
incompressible potential flow; the latter approximately 
characterizes the stagnation region even for compressible 
flows, except in the highly supersonic regime. 

Magnetic amplification beyond cb — 1, say, is not self- 
consistent, as our passive magnetization assumption fails. 
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Fig. 5.— Inverse polarity reversal layer (IPRL; defined by B z = 0; thick green curve) for substructure in the magnetic field far upstream 
of the sphere. Notations and flow are as in Figure 2. Here the magnetic field far upstream of the sphere is given by Equation (93), where 
9 0 = 7r/4, C = 0.25, oj = 6 tt, and t = 0. 


Instead, one expects a highly magnetized, self-regulating 
layer to form; its width <)' depends on the initial mag¬ 
netization level £b- For simplicity, consider an initially 
hot plasma, such that P th pu 2 for an inviscid flow, 
or P th pu 2 /Re for the Stokes flow. In this regime 
P th ~ P th , so according to Equation (97) 

s ^ Ub for negligible viscosity; 

I e ^ 2 for strong viscosity. 

For the compressible flows shown in Figure 6 , we use 
the approximate solutions of Keshet & Naor (2014, us¬ 
ing 7 = 5/3 and 0 = 0.52). In the supersonic regime, 
H is assumed to be upstream, so the jumps in B± and 
P th across the shock are taken into account. Here, the 
flow remains nonmagnetized until reaching the shock, at 
the so called standoff distance 6 = A; for high M 0 , the 
magnetized layer slightly shrinks as A diminishes. 

5.3. E Above Heliospheres and Magnetospheres 

The electric field in ideal MHD was shown to be a sim¬ 
ple function of the stream functions (see Equations (??, 
(46), (55)), so it is easy to infer E directly in various 
astronomical circumstances. As an example, in what fol¬ 
lows we consider the electric field outside a heliosphere 
(of some star) or magnetosphere (of some stellar ob¬ 
ject). We approximate the bow shock as weak (or ab¬ 
sent; for the possibility that there is no bow shock in 
front of the solar system, see McComas et al. 2012), and 
model the discontinuity (heliopause or magnetopause) 
using a Rankine half-body. These assumptions may be 
marginally invalid for the magnetopause of Earth, but 
roughly apply for the Sun’s heliosphere, and may be valid 
for other objects such as the induced magnetospheres of 
Venus and Mars. 

The Rankine half-body model was applied to the Sun’s 
heliopause by Parker (1961). Models for Earth’s mag¬ 
netopause are more sophisticated ( e.g Fairfield 1971; 
Tsyganenko 1989; Shue et al. 1997; Chao et al. 2002; Lin 
et al. 2010; Lu et al. 2011; Tsyganenko 2013, and refer¬ 
ences threin), and better calibrated by satellite observa¬ 
tions. However, no present analytic model captures the 
full complexity of the magnetopause structure, including 
the time dependence of the geodipole tilt direction and 


the stellar wind variability (Olson 1969; Zhou & Russell 
1997), so for simplicity we resort to Equation (99). 

The normalized potential of the Rankine half-body is 

<!> = -[ z + ^j , ( 99 ) 


where mixed coordinates are used, as in Section 4.1, and 
L m is the distance between the object (star, planet, or 
moon) and the nose of the discontinuity. This simplified 
potential flow is axisymmetric and incompressible. It ap¬ 
proximates the flow ahead of the discontinuity as homo¬ 
geneous and unidirectional, v = — uz , which is strictly 
valid only if the bow shock is weak or absent. The re¬ 
sulting discontinuity profile is approximately paraboloid 
at the nose, and quickly flattens to a g = 2 L m cylinder 
far downstream. 

The velocity field derived from Equation (99) is 


v = 



( QL 2 Jr 3 \ 

V—1 + ^/r 3 / 


( 100 ) 


The stream functions may thus be chosen as Ai = Rcj) 
and 


A2 = —- 
2 R L 


2 Li 


G-0 


= , ( 101 ) 

2 R 


where, as above, the length scale R is arbitrary. The 
thermal pressure is found from Bernoulli’s Equation (59) 

P,, = r.„ (x + ^ - ^=f) , (102) 


where P ram = pv 2 is the ram pressure, and x — 
Ph/Pram > 3/2 is necessary to guarantee a positive P th 
throughout the flow. 

The electric field outside the discontinuity is given, ac¬ 
cording to the axisymmetric flow Equation (55), by the 
simple analytic expression 


E 


Q_ 

QC 


1 v z B </> ^ 

\—v g B<p) 


(103) 
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Fig. 6 .— Head on (0 = 0) amplification of the magnetic pressure (left) and magnetic pressure fraction (right), plotted against the 
distance 5 from the nose of the unit sphere, for fields initially perpendicular to various steady flows: incompressible potential (solid blue; 
for P = pv 2 ), Stokes incompressible (dashed black; for P = 2 pv 2 /Re), and inviscid compressible (symbols with dotted lines to guide the 
eye) in the subsonic (M E {0.2,0.5,0.95}; green squares, small to large) and mildly supersonic (M E {1.2,1.4}; red circles, small to large) 
regimes. 


The magnetic field can be written in terms of elliptic 
functions (Roken et al. 2014). Close to the discontinuity 
g vanishes, so Equation (53) yields 



^ — BgVg/llS. 

B^g/L m I 

\-BgV Z /u) 


B Z + f Bg 
u 


( v e + 0(g)\ 

0 

\v z + 0(g)J 

(104) 


Near the discontinuity surface, where g is small, E is 
dominated by its g and z components, and B is domi¬ 
nated by the first term in Equation (104). In the cylindri¬ 
cal regime, v -+ (— uz ) and g -+ 2L m , so E# <C E e -C E z 
and « B 2 ~ 

The magnetic field diverges near the discontinuity. 
Hence, as in Section 5.2, we identify the region where 
€b > 1 as a highly-magnetized, self-limiting layer. The 
outer edge of the magnetized layer, defined, say, by 
es = 1 , can now be used to estimate the electromag¬ 
netic fields in the discontinuity region Cl, namely the he¬ 
liopause or the magnetopause. 

The width S of the highly magnetized layer varies 
slightly along the discontinuity, 


6 = 


5 + 3 cos 8 


1/2 


fL r 


(105) 


where 8 is the polar angle, taken with respect to the 2 
axis with the object at the origin, and we defined 
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The electromagnetic fields on Cl are thus given by 
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O(v^) ■ (108) 


As an example, consider the magnetopause of Earth. 
Satellite measurements (the 1999 Wind spacecraft tra¬ 
jectory; Wang et al. 2003) suggest that halfway between 
the bow shock and the plasma depletion layer, u — 
300 km s"\ B ~ 0.2 mG, B ~ (0.07, 0.17, -0.05) mG, 


and 



1 

2M\ 


— 0-1 , 


(109) 


X = 


Pr. 


1 

7 m 2 


~0.4, 


( 110 ) 


where P B = B 2 / (87 t) is the magnetic pressure, and Ma = 
v/(B 2 /Aitp) 1 / 2 is the Alfven Mach number. Hence the 
plasma is warm, and mildly magnetized with eg — 0.25. 

With these H parameters, the passive evolution of the 
fields in the Rankine half-body model Equations (105)- 
(108) yields (at 8 ~ 80°; Wang et al. 2003) S ~ 1.87?®, 
Bq — 0.35 mG, and Eq — 10 mV/m. These should 
be compared with the depletion layer thickness S dep — 
2.27?® (Wang et al. 2003) and with the magnetopause 
fields. Wind measured Bq = 0.35 mG (Wang et al. 2003), 
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but not the electric field; the CLUSTER spacecraft mea¬ 
sured Eq ~ lOmV/rn (but three years layer, at 9 ~ 70°; 
Panov et al. 2006). The good agreement between mea¬ 
surements and model is largely coincidental, as the latter 
is crude, pertains only to the draped component, and fails 
to account for the complexity of the magnetosphere, the 
variability of the solar wind, and physical effects such as 
magnetic reconnection. 

6 . SUMMARY AND DISCUSSION 

We have shown that the stream functions A or path 
functions A, often used to describe non-diffusive flows, 
provide a simple, local (in A or A space), and general so¬ 
lution to the problem of evolving electromagnetic fields in 
a highly conductive plasma, or equivalently propagating 
a length element l, in an arbitrary flow. Explicit expres¬ 
sions were provided, in both steady and time-dependent 
flow regimes, for the evolving magnetic field B (Equa¬ 
tions (29) and (45)), electric field E (Equations (27) and 
(46)), and l (Equations (28) and (44)). 

The analysis, illustrated in Figure 1 for a steady flow, 
uses A or A as coordinates. Vectors are thus projected 
onto a reciprocal basis (Equations (14) or (44)), in which 
the vector components are conserved along streamlines 
or pathlines. It is in this basis that the analysis be¬ 
comes manifestly local, and thus conceptually and com¬ 
putationally simpler than previous solutions of the ideal 
MHD equations. The results are more general than avail¬ 
able until now because the closed-form expressions for B 
and E can be applied ”as is”, given A or A, for arbitrary 
magnetic field boundary conditions specified on any hy¬ 
persurface H. The analysis is valid for an arbitrary ideal 
MHD flow. 

For steady flows, a third variable q is needed, in ad¬ 
dition to the two stream functions, in order to map the 
volume spanned by the flow. Our analysis utilizes the 
freedom in the choice of this q\ see Equations (22) and 
(25). Evolving B becomes considerably easier for certain 
choices of q 1 in particular if the travel time parameter At 
is measured or easily computed; other choices require, in 
general, one integration (Equation (23)). Conveniently, 
the evolution of E is entirely fixed by the stream func¬ 
tions alone, given boundary conditions on H. 

For time-dependent flows, the analysis is conceptually 
simpler because the three path functions entirely deter¬ 
mine the flow. The freedom in the choice of a fourth 
variable Q , needed to fully map spacetime, is less im¬ 
portant. Here, no integration is necessary when A a are 
known. The analysis is demonstrated for an arbitrary 
ID flow (Equations (86) and (87)). 

For both steady and time-dependent flows, computing 
the electric field is particularly simple in our method, as 
it only requires a differentiation of the stream or path 
functions. This is much simpler than first computing the 
magnetic field, and then using Ohm’s law to derive E. 
because finding B typically involves a numerical integra¬ 
tion (unless the path function are known). 

The general expressions for the electromagnetic field 
allow us to quickly reproduce and generalize known mag¬ 
netization solutions, such as for a potential incompress¬ 
ible flow around a sphere (Equations (68) and Figure 2), 
and to derive new solutions for magnetization in vari¬ 
ous flows, such as an arbitrary axially symmetric flow 
(Equation (53)), the incompressible Stokes flow around 


a sphere (Equations (79) and Figure 3), and the com¬ 
pressible flow in front of a blunt axisynnnetric object 
(Equation (96) and Figure 6). 

The magnetization in front of an axisynnnetric blunt 
object (Equation (56)) is of special interest in various as¬ 
tronomical systems, such as planetary bodies moving in 
a stellar wind, a stellar system moving through the ISM, 
or clumps and bubbles moving through the IGM. The 
analysis simplifies considerably if the flow along the axis 
is approximately irrotational (as C in Equation (57) van¬ 
ishes) . The simple relation thus obtained between B and 
v can be used to directly infer one quantity, if the other 
is known. Here, the magnetic pressure (Equation (96)) 
and magnetic pressure fraction (Equation (95)) increase 
toward the stagnation point, roughly as i> _1 (see Equa¬ 
tion (97)), as illustrated for various flows in Figure 6. 

Viscosity and compressibility are seen (in Figure 6) 
to enhance the magnetization, leading to a thicker mag¬ 
netized layer. In contrast, the magnetic pressure frac¬ 
tion appears to be quite insensitive to compressibility. 
Crudely, the thickness of the magnetic layer depends on 
the upstream magnetization is through Equation (98). 

One advantage of our approach is that the results 
are directly generalized for nonhomogeneous electromag¬ 
netic field configurations upstream. Indeed, given some 
spacetime-dependent Bi„ on H , the analysis remains in¬ 
tact, provided that B is determined by Equation (88). 
This is illustrated in Figures 4 and 5 for Fourier compo¬ 
nents of Bi n in the stellar wind in front of a planetary 
body. As a concrete example, we derive a semi-analytic 
expression (Equation (92)) for the position of the IPRL 
behind a planetary body, and demonstrate (in Figure 5) 
the substructure it acquires due to upstream magnetic 
variability. 

Another advantage is that the electric field can be 
determined analytically in various flows. As an illus¬ 
tration, we compute E (Equation (103)) for a Rankine 
half-body model, approximately applicable above some 
heliospheres and magnetospheres. Modeling the tangen¬ 
tial discontinuity as a self-limiting magnetized layer, we 
estimate the draped electromagnetic component of the 
corresponding heliopause and magnetopause (Equations 
(105)-(108)). 

Our approach directly links the properties of electro¬ 
magnetic fields to the properties of path functions and 
stream functions. This connection, unexplored previ¬ 
ously as far as we know, may thus help shed light on 
the general properties of both electromagnetic fields and 
path or stream functions, in complicated flow or field 
configurations. Even if only one stream function is (or 
only 1 -2 path functions are) known, our approach may 
still provide valuable constraints on the electromagnetic 
field evolution. Vice versa, if the electromagnetic fields 
are known, the stream (or path) functions may be con¬ 
strained, or even determined. 

An important advantage of our analysis is its direct ap¬ 
plicability to hydrodynamic numerical simulations that 
are already based on stream functions or path functions. 
Indeed, our results (in particular Equations (25),(27) or 
Equations (45),(46)) can be immediately incorporated 
into such simulations, in order to follow the passive evo¬ 
lution of electromagnetic fields. Moreover, our results 
suggest that the back-reaction of electromagnetic fields 
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on the flow could be easily added to such a code, thus 
leading to an efficient MHD simulation. 

Consider, for example, a putative path function simu¬ 
lation for a 3D time-dependent flow. It is most natural to 
choose the hypersurface H as the 3D volume at some ini¬ 
tial time; this corresponds to imposing initial conditions 
on the magnetic held, in addition to the initial condi¬ 
tions on the how. As the grid coordinates correspond 
to the path functions, the 3D magnetic configuration is 
simply frozen onto the grid. Thus, the magnetic held 
does not need to be evolved, and its back reaction on the 
how can be easily computed. Implementation schemes 
are deferred to future work. 

Hydrodynamic simulations based on stream functions 
or path functions have been slowly evolving over the past 
few decades. Stream function codes have been developed 
for steady hows in 3D ( e.g Pearson 1981; Clrungyuan 
& Dulikravich 1986; Sherif & Hafez 1988; Beale 1993; 
Knight & Mallinson 1996; Elshabka & Chung 1999). 
Path function simulations have been so far presented for 
time-dependent hows in ID (Hui 2007) and in 2D (Loh 
& Hui 2000). 

Like any Lagrangian code, such simulations have the 


advantage of precisely conserving the mass flux, match¬ 
ing the resolution to the evolving how, and sharply re¬ 
solving discontinuities. Employing path functions has 
the additional advantage of converting the hydrodynamic 
equations into a precise, geometric, conservative form. 
This bears considerable numerical advantages, and guar¬ 
antees the correct deformation of each grid cell, although 
strong cell distortion e.g., in vortical how, must still be 
addressed. This has been explicitly shown in 2D (Loh 
& Hui 2000), and partly developed for 3D (Van Roes- 
sel & Hui 1991). An MHD simulation based on path 
function will thus have important advantages in both 
hydrodynamic and electromagnetic sectors, and would 
supplement present MHD codes. 
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APPENDIX 

A. THE RELATION B <x pi PRESERVES THE MHD EQUATIONS 

We evolve the electromagnetic fields using Equation (8), or equivalently the Helmholtz equation (6), so it is useful 
to prove that the resulting B satisfies Gauss’ law (2) and the convection equation (4). We assume that the Helmholtz 
equation and the continuity equation (3) are satisfied everywhere, and that Gauss’ law holds on H. 

First note that 

— - V x (« x B) = — - (B • V) v + (v • V) B + B (V • v) - v (V • B) = (-«) (V • B) , (Al) 

ot at 

where we used the Helmholtz and continuity equations in the last step. This ensures that if Gauss’ law is satisfied, 
then the convection equation (4) is guaranteed, and vice versa. Taking the divergence of this equation yields 

diy dt B) =-(y-v)(y-B) . (A2) 

With the initial conditions V • B = 0 assumed on H , the trivial solution Equation (2) spans the entire flow connected 
to H. Gauss law is thus satisfied, and, according to Equation (Al), so is the convection equation. 


B. THE INTEGRAL {d\M) r FOR A STEADY INCOMPRESSIBLE POTENTIAL FLOW AROUND A SPHERE 

The integral in Equation (69) can be solved by expanding the integrand about A = 0. This expansion is useful for 
analyzing the flow in the close vicinity of the axis of symmetry or of the sphere, where A is small. The result is 
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r is the Euler Gamma function, and 2 -F 1 is the hypergeometric function. In particular, 
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